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1. Introduction 

In QCD, for many observables of interest, the only calculational approach that is 
currently available is a perturbative expansion in powers of the strong coupling as- 
Quite often the coefficients of this expansion are well behaved and reasonable accu¬ 
racy can be obtained given the first two or three terms of the series. 
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However there exist whole classes of observables for which the coefficients grow 
rapidly order by order. Typically this happens in problems involving two or more 
scales, in which case the coefficient of o;” may involve up to 2n factors of the loga¬ 
rithm, L, of the ratio of scales. These logarithms are generally associated with soft 
and collinear regions of phase space. When L is large it compensates the smallness 
of ttg and reliable predictions can be obtained only by including an inhnite number 
of terms. While signihcant technical difficulties prevent an exact determination of 
high-order terms, it is sometimes possible to calculate the coefficient of their dom¬ 
inant (and sub-dominant) logarithmic parts, order by order, and resum this subset 
of terms. In many cases this is sufficient to obtain a reliable prediction. 

One class of observables that has been extensively studied in this context is that 
of event shapes [1-12] and jet rates [13-15]. These measure geometrical properties of 
the hadronic energy flow in the hnal state of a given reaction, e.g. in e+e“ —hadrons, 
and are widely used in tests of QCD, the measurement of ag, and studies of the 
interface between perturbative and non-perturbative physics. When the value V of 
the observable is small, the perturbative expansion of its distribution is dominated 
by terms (o;” In^""”^ 1^)/!^ a resummation is necessary. 

The determination of the order ag and aj contributions to such observables is 
usually fairly straightforward. One writes a computer subroutine which calculates the 
observable for an arbitrary set of 4-momenta and links it with a program which codes 
the hxed-order real and virtual matrix elements for a given process. The ‘hard’ part 
of the problem is to determine the matrix elements and code them into the program 
so as to obtain a numerically efficient cancellation of real and virtual contributions. 
But this part of the problem is common to all infrared and collinear safe observables 
and is ‘ready-solved’ in the form of a range of publicly-available programs, such as 
EVENT2 and DISENT [16], and also many others, see for example [17-22]. 

On the other hand the determination of resummed distributions of event shape 
observables has up to now proved more labour-intensive. This is despite the fact 
that many aspects of the calculation are actually fairly similar from one observable 
to another — for example one can generally make the same set of approximations 
concerning the multi-particle matrix-element, writing it as a product of matrix ele¬ 
ments for independent emissions. However it is also necessary to relate the value of 
the observable to some condition on the ensemble of emissions. For a resummation 
to be feasible analytically, one must express this condition as a product of factors 
each involving a single emission. It is this part that requires a separate analytical 
calculation for each observable. 

To explain in more detail what we mean, it is useful to consider the thrust T in 
e’''e“ (dehned in appendix B).^ In the two-jet region its value is given by a sum over 

^The thrust is slightly unusual in that the two-jet region (where the resummation is needed) 
corresponds to T ^ 1, but this is just a matter of definition. 
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contributions from individual emissions: 


T = (1.1) 

* ^ 

In performing the resummation one determines the cross section for 1 — T to be 
smaller than some value, say 1 — T < r. Using an integral transform this condition 
can be written as 

e(r-(i-r)) = |. (12) 

which is a product of factors for individual emissions, just as required. 

The thrust is a particularly simple example, because it is easy both to determine 
its dependence on multiple emissions, eq. (1.1), and to factorise this dependence, 
eq. (1.2). Other variables are less trivial. For example for the jet-broadenings the 
determination of the dependence on multiple emissions involves subtle issues of hard- 
parton recoil, and the factorisation requires at least two integral transforms [4]. Some 
multi-jet event shapes, such as the thrust minor in the three-jet limit involve as many 
as hve additional Fourier transforms [8]! This makes the resummation of event shapes 
quite tedious (and sometimes even error prone). Furthermore for some variables, for 
example the thrust major and oblateness in it is not even possible to write down 
a closed expression for the dependence of the observable on a general multi-particle 
conhguration, rendering an analytical resummation unfeasible. 

So in this article we take a step towards the development of a systematic semi- 
numerical approach for calculating resummed distributions for a large class of observ¬ 
ables, essentially those for which double logarithms exponentiate. The probability, 
S(n), of a suitable observable having a value smaller than v is given by an expression 
of the form 


J:{v) = exp[Lgi{asL) + g 2 {asL) ^ -], L = \nl/v. (1.3) 

The terms in Lgi{asL) are known as leading, or double logarithms (LL, DL), and 
start at 0{asL‘^). The terms in g 2 {asL) are known as next-to-leading, or single 
logarithms (NLL, SL), and start at 0{asL). The double logarithms arise from 
a veto of a soft and collinear portion of phase space and are straightforward to 
calculate. The single logarithms arise from a variety of sources. For global variables 
[12] most of the sources {e.g. from the running coupling, hard collinear emission, 
angles between hard partons) can be understood by considering the phase space 
for just a single emission. What remain are single logarithms associated with the 
observable’s dependence on multiple emissions. It is this contribution that we address 
here. 

Our method involves hrst considering a ‘simple’ reference variable, in which the 
dependence on multiple emissions is trivial, leading to a straightforward analytical 
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resummation. One then numerically relates the resummation of the simple variable 
to that of a more complicated variable of our choosing. For this to work the two 
variables must have the same structure of double logs. 

So in section 2 we show how to relate a pair of variables with the same double 
logs, and give a procedure for designing a suitable ‘simple’ reference variable. We 
then show how this can be implemented in a Monte Carlo algorithm for calculating 
the single-logarithmic function relating the resummation of the two variables. 

In section 3 we show that this method reproduces known results, with examples 
from e+e“ —2 jets and e+e“ —3 jets. 

Finally in section 4 we present new results for some widely studied variables in 
e’''e“ —2 jets: for the thrust major and the oblateness for which no resummation 
had ever been performed up to now (and which are probably beyond analytical 
treatment) and for the Durham three-jet resolution parameter for which resummed 
results existed but not to full NLL order [13,14]. 

2. General procedure 

This section sets out the principal elements of our approach. It is divided into several 
parts. First we see how, at NLL accuracy, to relate the distributions of two variables, 
starting from knowledge of the probability of distribution of one variable given the 
value of the other. 

We then see how to design a ‘simple’ reference variable for which the analytical 
resummation is trivial and examine how to calculate the probability distribution of 
a more complicated variable given the value of this simple one. 

The section finishes with the consideration of the subtle but practically important 
issues of how to ensure that our hnal numerical calculation is free of contamination 
from next-to-next-to-leading logarithmic (NNLL) corrections, and how to determine 
its expansion to second order in ag (needed for the subtraction of doubly counted 
terms when matching to fixed order results). 

2.1 Relating the resummations of two observables 

Suppose we have a ‘complicated’ observable V and a ‘simple’ observable W with 
the same all-orders double logarithmic structure (terms in the exponent). 

We introduce the probability distribution P(u|us) for the value v of the complicated 
observable, given a value Vg for the simple one. 

Writing the distribution of the simple observable as Ds{vs)/vs-, one obtains the 
distribution of the more complex observable D{v)/v through the following convolu¬ 
tion: 

D{v) = v [ — Ds{vs) P{v\vg) . (2.1) 

J Vg 
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If P{v\vs) is dominated by the region v ^ Vg (this will follow naturally from the two 
variables having the same double-log structure) then one may expand the distribution 
Ds{vs) using eq. (1.3), 


Ds{Vs) 


Ds{v) exp 


-R' In - + C> 

Vg 



( 2 . 2 ) 


where R', which is related to the (differential) phase space for emissions, is dehned 
as 


R' 


dlnS 

dL 


(2.3) 


Neglecting NNLL terms, ttgL” this gives 


D{v) = vDg{v) 


Vs 


(2.4) 


As we shall see shortly, for suitable observables (as dehned in appendix A), to NLL 
accuracy vP{v\vs) depends only on the ratio v/vg and on the phase space for emissions 
(through R'): 


= + (2.5) 

where we have introduced a rescaled probability p{x, R') for the ratio v/vg to be to 
equal to x (with measure dx/x). As long as v and Vg are of the same order, the 
ambiguity of whether R' should be evaluated at v or Vg corresponds to a subleading 
effect since R'{v) = R'ivg) + O hiv/vg). This allows us to write 


D{v) = Dg{v)P{R!), 


( 2 . 6 ) 


where 

R{R') = [ ^e-^'^^^p{x,R') . (2.7) 

J X 

A relation analogous to eq. (2.6) holds also for the integrated distributions 

S(u) = J:g{v)P{R'), (2.8) 

as can be seen by differentiating and noting that extra term relative to (2.6) is 
subleading: 

TE(„)=(^(fl') + |^)TE». (2.9) 
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2.2 The ‘simple’ observable 

So far the simple observable has been left as a fairly vague concept. Here we outline 
a concrete procedure for constructing it. 

Given a Born conhguration consisting of hard momenta {pBom} together with 
an arbitrary set of soft and collinear emissions /ci,... , /c^, a general way of dehning 
the simple observable is as follows: 


^({PBorn}) ^1) ^2) ■ ■ ■ t ^n) 

max[V{{pBorn},kl),V{{pBorn},k2), . . . , H({pBorn},/^n)] • (2.10) 


This dehnition is unsafe with respect to secondary collinear branching of the soft and 
collinear emissions. However for the ‘suitable’ class of event shapes under considera¬ 
tion, secondary collinear branching is irrelevant except insofar as it is responsible for 
determining the scale of the coupling [1] — accordingly to calculate the resummation 
of the simple observable and we consider only independent emissions from the 

Born conhguration and account for collinear splitting by directly setting the scale of 
the coupling.^ 

Having dealt with this issue, the resummation for Vg is then simple, because the 
condition for the value of the observable to be smaller than Vg factorises straightfor¬ 
wardly, 

n 

Q{Vg - Vg{ku k2,..., kn)) = n - y{k;)) , (2.11) 

i=l 

where for brevity we have dropped the Born momenta from the arguments of V. 

To illustrate this we shall consider the class of global observables in 
2 jets. Emissions can be taken to be independent, as shown in [2,4], giving the 
following expression for their probability distribution: 


dP(ki,... ,kn) 


n! ki 271 tt 

1=1 


V hi 4 



( 2 . 12 ) 


where hi, rji and (f)i are respectively the transverse momentum, the rapidity and the 
azimuthal angle of h with respect to the qq axis. Virtual corrections are given by 
the factor with 


Re 


f ^^dhj d(j)CFC(hkt) ^ 

J hi 277 77 




(2.13) 


where the e indicates that some regularisation must be applied to both the real 
and virtual parts. In the above formulae, the coupling is dehned in the gluon 

^One could perhaps devise a collinear-safe definition which for our purposes is equivalent to the 
above one. This might involve the application of a clustering algorithm to emissions before using 
eq. (2.10). 
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Bremsstrahlung scheme [23] and the —3/4 in the limit on the rapidity accounts for 
the hard part of the Pqq collinear splitting function (see for example [4]). Combining 
these expressions with eq. (2.11) gives 

Es(t;J = exp (-i?^(nj) , (2.14) 


with 


Rs{Vs) — Cp 


dkf d(j) as{kt 




kf 271 


71 




dr]Q{V{k) — Vs 


(2.15) 


In processes with three hard partons the resummation of the simple variable 
also proves straightforward, with Rg just involving a sum over dipoles [8-11]. On the 
other hand, with 4 or more hard partons subtleties will arise from the large angle 
region as can be seen from [24].^ 


2.3 The probability P{v\vs) 

One of the advantages of the above choice for the ‘simple’ observable is that it is 
easy to design a Monte Carlo algorithm to construct, for a given R', an ensemble of 
conhgurations which all have the same value of Vg and known (equal) weights. 

Fixing Vg implies taking only those conhgurations where one of the emissions, 
j, satishes V{kj) = Vg, and all others i ^ j satisfy V{ki) < Vg. Using eq. (2.12), 
renumbering emissions such that j is always 1, and dividing by the total probability 
of the simple variable having value Vg, we obtain the following expression for the 
probability distribution of conhgurations given Vg-. 


V{h)) 

where 

rl R 

R'givg) = -Vs^ = R'ivg) + O . 

dVs 

It is now straightforward to construct P(n|n 5 ): 


-R,+Rg{vg) « 

dP{ki,... ,kn\vg) = " ^ X (5(ns - U(A;i)) TT0(n,- 

Rsi'Vs) (n-1)! 


2=1 


kti 27r 


7T 


kti 


(2.16) 


(2.17) 


P(n|nJ 


oo „ 

/ dP{ki,... ,kn\vg)6{V{ki,... ,kn) - v). 

n=l 


(2.18) 


^Though it remains to be seen whether the extensive technology developed in [24] accounts for 
the non-global effects that are to be expected for the observables that are considered there. 
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2.4 Eliminating subleading effects 

As they stand, the above expressions contain not just single-log terms but also many 
(sometimes spurious) subleading contributions. For example while in eq. (2.15) the 
‘3/4’ is necessary in order to obtain the correct single (NL) logs, in eq. (2.16) it leads 
only to NNLL terms. These subleading terms are awkward for a variety of reasons. 
Firstly there is currently no way of guaranteeing their correctness. Secondly when 
matching to hxed-order predictions they lead to some technical difficulties because 
one needs a good knowledge of their expansion to 0{ag). Finally one notes that 
the expressions of the previous section contain integrals over as{kt) down into the 
infrared. The Landau pole in Ug (or whatever other structure occurs in higher orders) 
will then introduce a power correction ambiguity in one’s answer, requiring that one 
make an arbitrary non-perturbative ‘choice’ concerning its treatment. Accordingly 
it has become standard to give predictions containing nothing but leading and NL 
logs (though we note the alternative philosophy advocated in [25]). 

Truncating the expressions for Rg and R' at NLL order is easy, since they are 
computed analytically. However ensuring that T is purely a function of R!^ required 
in order for it to have just NLL terms, turns out to be the most subtle part of our 
semi-numerical resummation approach. We present two ways of doing it. 

Method 1 (general). One way of eliminating subleading effects is by numerically 
taking the limit of Og 0, while keeping Oglnl/n hxed (he. hxed R!\ While 
conceptually straightforward and very general, such an approach turns out to be 
numerically awkward: the NNLL contamination in one’s answer is of order Og, and for 
hxed /?', In 1/n scales roughly as l/os- So if for argument’s sake we want to eliminate 
subleading effects to within a percent, then we need to consider values of v of the order 
of (since we are dealing with orders of magnitude of the logarithm, we ignore 

the difference between and 10“^°°). Unfortunately standard algorithms for 

calculating event shapes generally work only for n 3> where e is the relative hoating 
precision (usually 10“^® in double precision) and p is some algorithm-dependent 
power of order 1. 

So in order for this approach to work we need to hnd an algorithm for calculat¬ 
ing the observable that is free of sensitivity to the limited accuracy of the hoating- 
point arithmetic.^ Additionally the algorithm should work with a representation of 
4-momenta that does not suffer from the limit on the smallest representable floating¬ 
point number (in standard double precision the smallest number that can be rep¬ 
resented is of the order of 10“^°®). Both these issues restrict the generality of the 
approach, and imply a certain amount of analysis of one’s observable so as to un¬ 
derstand the elements fundamental to its determination in terms of a set of soft and 

^An alternative approach, whose feasibility we have not investigated in detail, would be to use 
arbitrary precision arithmetic. This would entail a significant penalty in terms of computing time 
needed, but might nevertheless be an approach worthy of investigation. 



collinear momenta. This will be illustrated in detail in section 4.3 when we come to 
discuss the 2-jet rate for the Durham algorithm. 


Method 2 (event-shape specific). Fortunately there exists a much simpler ap¬ 
proach which works in almost all cases. Event shapes (but not jet rates) have the fol¬ 
lowing useful property: given an ensemble of soft and collinear momenta {fci,... , kn\, 
if we vary the {ki\ such that the {V{ki)} are kept constant and the {ki} all remain 
collinear to the legs to which they were originally collinear, then V{{ki}) also remains 
constant. In eq. (2.18) this allows us to carry out the integration over rapidity ana¬ 
lytically and replace it with a sum over the different hard legs to which an emission 
may be collinear (in calculating JF the large-angle region gives only NNLL terms). 

So we divide the phase space R' into contributions R'^ coming from N different 
hard legs 


N 

«'(«) = E F(0 . (2.19) 

e=i 

and write 


dP(ki, • • • , kn\vs) 


y y — Re-\-Rs{Vs) 

X ^(ti,-t'i), (2.20) 

i=2 




i=l 


Vj 271 


where the use of R rather than R denotes the fact that emissions are characterised 
by the leg ii to which they are collinear, their azimuthal angle 0* and the value of Uj. 
The transverse momentum and rapidity are no longer individually dehned — they 
have been integrated over while keeping vi constant. 

To eliminate subleading logarithms, we now exploit the infrared safety of event 
shapes, which implies that the emissions that contribute significantly to V{{ki}) are 
those for which Uj ~ ui. Accordingly we can throw away logs of vi/vi (both in the 
real and virtual parts) to get 


dP(ki,, kn\vs) = 


V. e 


_D/ rvs 

J V 


R' (n-l)! 


n 

2=1 


dvi d(j)^ , 
Vi Itt 


ti 


8{vs 


V\ 


ne(. 


( 2 . 21 ) 


i=2 


where the are all evaluated at Vg. By truncating their analytical expressions at SL 
level, we obtain an answer for dP^... |us) that has only single logarithmic dependence 
on Vg. Finally to get the equivalent of eq. (2.18) we have to remember not only to 
integrate over the 0* and Vi but also to sum over all the ii. 


p(-,R)=vf:[ 

dP{ki,... ,kri[ 


£l,... , In 


, K) - v) 


( 2 . 22 ) 
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We have written the answer for p{v/vs, R') rather than for P(n|ns) so as to emphasise 
that we are free to choose Vg as we like, and that the result is a function only oi v/vg 
and R', free of subleading logarithms (we recall that at NLL order, R' itself depends 
only on a^L). In practice there is still some contamination from power suppressed 
contributions associated with the observable itself, terms of relative order with q 
some observable dependent power of order 1. However it is usually perfectly feasible 
to choose values of Vg such that 1, while keeping the relative rounding errors 

associated with the numerical calculation of the observable, s^/vg, also negligible. 

For the calculations whose results are given below, the following Monte Carlo 
algorithm has been used to generate conhgurations according to the probability dis¬ 
tribution in eq. (2.21). One starts with i = 1 and yi = 1, then: 

1. One chooses a leg ii randomly, the probability of leg i being RyR'. 

2. One chooses uniformly a random azimuthal angle, <pi, with respect to the leg. 

3. One then chooses the transverse momentum ku and rapidity rji with respect 
to the leg, such that the emission is soft, collinear and satisfies V{ki) = yiVg. 
Within these constraints one has total freedom in one’s choice of transverse 
momentum and rapidity. 

4. One chooses In^/j+i such that |/j+i < yi, with a random distribution propor¬ 
tional to {yi^i/yi)^' reflecting the fact that the phase space to produce an 
emission with yvg < Vk < {y + Sy)vg is R'5y/y. 

5. If i/j+i < e 1, with e some arbitrary small cutoff, one stops. Otherwise one 
replaces i ^ i + 1 and goes back to step 1. 

2.5 The coefficient of R'^ 

For the purposes of matching with fixed order calculations it is necessary to have the 
expansion of R{R') in powers of R' 

OO 

R{R') = l + ^RnR"" (2.23) 

n=l 

(from which one can deduce the expansion in terms of powers of Os-^)- Since matching 
is currently usually carried out to NL order one needs the first two terms. As a 
consequence of our definition of the ‘simple’ observable, the first term is automatically 
zero. The second term can be obtained by expanding eqs. (2.7) and (2.22) to give 

V{ki) = Vg . 

(2.24) 


^2 = - 


Y. 


lim 

asL^O 


R'^ 


d(l)\d(j)2 

(27r)2 


V2 Vg 
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We specify the limit -R' —0 for the ratio of R'^/R' because for some more complicated 
observables, such as Tze in DIS [5], RyR' is not a pure number but rather itself a 
function of R'. We also note that despite the appearance of Vg in the integral, is 
of course independent of Vg. 

Finally we point out that if necessary y.g. for NNLO matching) one can in 
a similar way derive corresponding expressions for yet higher order terms in the 
expansion of T. 


3. Reproducing known results 


An important test of our method is that it reproduces known analytical results. This 
section shows results for a variety of previously calculated observables in conhgura- 
tions with 2 and 3 hard jets. For dehnitions of the observables, the reader is referred 
to appendix B. 


3.1 e+e —^ 2 jets: thrust and broadenings 

Here we recover the known results for the thrust T and the jet broadenings (total 
Be and wide B^). In each case we have tested both methods for the elimination of 
subleading logs (section 2.4) and verihed that they give identical results. 

For the thrust distribution the ‘simple’ observable is given by 

Tg = max I I , (3.1) 


The lowest line of figure 1 shows the comparison between the function Rr and the 
corresponding exact result [1], 


p-isR' 

^ “ F(1 + R) ■ 

In the broadening case we start from 


Bg = max 

i 


I kti 

~Q 


(3.2) 


(3.3) 


The Monte Carlo procedure gives Bl and Br (see (B.5)) and then Bt and Bw (see 
(B.6)). The comparison between numerical and analytical resummation [4], 


RbAR') = 




q—'yeR' 


dx /1 + x\ ^ 


F(l + i?')Ui 

p-1eR' 


F2(l + R'/2) 
is given by the top two lines of hgure 1. 


dx f 1 + x\ ^ 


x^ 


(3.4) 
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R' 

Figure 1: T{R') for the thrust and the two broadenings. The lines are the analytical 
results, while the points are the numerical results. 

3.2 e+e“ —^ 3 jets: thrust minor (a.k.a. X^out) 

The procedure can also be applied to multi-jet event shapes, such as and D 
parameter, which have been studied recently in the near-to-planar three-jet region [8, 
9]. In particular, we would like to compare the numerical results to the analytical 
resummation of distribution, which is much more involved than that of the D 
parameter due to hard parton recoil effects. 

In the near to planar three-jet region T ~ Tm S> T^, a three-jet event consists 
of a hard quark-antiquark-gluon system accompanied by soft secondary partons. We 
denote the hard partons by pi,p 2 and ps with p? > p® > Ta ^ 2, 3 the 

configuration in which the gluon momentum is ps. As discussed in [8], the function 
depends on the colour configuration of the hard underlying system. 

The simple observable Tm,s is determined by (2.10). In terms of soft parton 
momenta this gives Tjn^s = where ki^ is the out-of-event-plane 

momentum component of emitted parton i and, due to the recoil kinematics, W is 
4 or 2 respectively according to whether or not emission i is in the same hemisphere 
as the most energetic hard parton (as was shown in [8] and can easily be determined 
numerically as well as analytically). 

The analytical result for [R') has been computed in [8].® We do not reproduce 
®Note that in [8] T represented a different quantity. 
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R' 

Figure 2: for the thrust minor distributions in the three different colour con¬ 

figurations. The lines indicate the analytical results, while the points are the numerical 
results. 

here its explicit form, since its various components involve up to as many as five 
nested integrals, but we plot in figure 2 the analytical results as a function of R! 
together with obtained using the numerical procedure. 

4. New results 

We now exploit our method to compute the function JF for some observables for 
which an analytical expression for the resummed PT distribution has not so far been 
found: the thrust major, the oblateness and the Durham three-jet resolution. In 
the hrst two cases we suspect that it may not even be possible to hnd analytical 
expressions. 

4.1 Thrust major 

The thrust major, T^, is dehned in (B.2). In order to compute Ttm have as the 
‘simple’ observable 

Tm,s = 2 ■ max . (4.1) 

Our numerical procedure gives the function shown in hgure 3. The re¬ 

summed Tm distribution is then given by (2.8), where Ss(Tm), dehned generally in 
(2.14), is 
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Figure 3: The function J-{R') for the thrust major, the oblateness and the Durham 
three-jet resolution. 




Rs(Tm) — ‘2Cf 



dk^ asih) 
kt tt 




(4.2) 


i.e. it is expressed in terms of a radiator, Rs(Tm), which is identical to that for Bt 
and B]y. Its analytical expression to SL accuracy has been already computed in [2] 
and is recalled in Appendix C.l. 

We then match the above resummed result (using the log i?-matching scheme [1]) 
with the exact O (a|) results computed with EVENT2 [16], to give the curve shown 
in hgure 4. 


4.2 Oblateness 

The ‘simple’ observable (4.1) can be also exploited to compute B{B!) for the oblate¬ 
ness, defined as the difference between the thrust major and the thrust minor (see 
(B.4)). 

The function B{R!) for the oblateness (also shown in figure 3) behaves quite 
differently from that for the other variables considered so far, since it increases rather 
than decreases with increasing R'. 

There is a simple reason for this difference. For most variables, adding an extra 
emission to the ensemble can only increase the value of the observable, largely because 
most observables essentially involve a sum of positive dehnite quantities. 
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Figure 4: The full SL-resummed distributions for the thrust major and the oblateness 
including second order matching (log i?-scheme), shown for t/s = Mz- 

The oblateness is unusual in that it is a difference of two quantities. For a single 
emission the thrust minor is zero and so the oblateness is equal to the thrust major. 
But as soon as one considers a second emission, there exist configurations for which 
the thrust minor is almost equal to the major, leading to a value for the oblateness 
which is smaller than that of the simple observable. 

This is illustrated in figure 5 which shows the probability distribution for the 
ratio of the true variable to the simple one, for both the thrust major and the 
oblateness. The thrust major exhibits a sharp cutoff for v/vg < 1. The oblateness 
on the other hand has a roughly power-like tail extending to arbitrarily small values 
of v/vg, originating from events with nearly identical major and minor projections. 
Since eq. (2.7) involves the ratio {v/vg)~^\ if a sufficient fraction of events has v < Vg, 
then T{R') can be larger than 1. 

However T being larger than one is not the only consequence of the tail at 
small vjvg. if B! is too large then the factor {y/vg)~^' can completely compensate 
the smallness of P{v\vg), and the integral (2.4) then diverges. The start of this 
divergence is clearly visible for the oblateness in hgure 3. We shall refer to the 
position of the divergence as 

In the case of the oblateness, simple considerations based on two-gluon conhg- 
urations suggest that = 2, corresponding to a tail of vP{v\vg) proportional to 
(v/vg)^. This does not quite correspond to what is seen in hgure 5, and it is not clear 
whether the actual vP{v\vg) has a different power for the tail, or simply some extra 
logarithmic enhancement. 
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Figure 5: Distribution of the value of v/vg for the oblateness and the thrust major. 

We emphasise that the divergence in T is not a particularity of our approach. 
Recently for example, such a problem has been extensively discussed in the context 
of the jet broadening in DIS [6], and it is in general present for all observables which 
can have zero value even in the presence of emissions. Physically, what happens is 
that for i? ~ there is a transition from the normal double-log Sudakov suppression 
to some other behaviour {e.g. a suppression proportional to n^). Such a transition 
is beyond representation in terms of a pure SL function, leading to a breakdown in 
the hierarchy of leading and subleading logs. This is seen for example in the fact 
that NNLL terms are even more strongly divergent that the NLL ones. In analytical 
resummation approaches, with considerable extra work, it is sometimes possible to 
include an appropriate subset of subleading terms so as to bring the answer back 
under control over the whole range of the R\ 

On the other hand as long — R' is of order one (formally 3> it can be 

shown that the divergence can be ignored [6]. Therefore, so as to be sure of remaining 
in the region under control, we suggest that the oblateness distribution be studied 
only for R'{0) < 1. This is illustrated in hgure 4, which shows the full resummed 
matched distribution for the oblateness (based on an expression analogous to (4.2)). 
The point R' = 1 is just to the right of the peak region — accordingly, it should be 
possible to compare most, though not all the data with the resummed predictions. 

4.3 Durham three-jet resolution 

Now we wish to calculate the multiple emission single-logs for the distribution of 
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the three-jet resolution parameter, 1 / 3 , in the Durham algorithm (the corresponding 
integrated distribution is also referred to as the two-jet rate). The leading logs, and 
some subleading logs (those necessary for NLL accuracy in the exponentiated answer 
rather than the exponent) were calculated in [13]. The NL logs associated with the 
appropriate scheme choice for as were given in [14]. However those coming from the 
non-trivial dependence on multiple emissions have yet to be calculated. 

Jet rates differ from most event shapes in that they do not satisfy item 3 of the 
conditions for applicability given in appendix A, which was required for a straight¬ 
forward elimination of subleading logarithms, using the event-shape specihc method 
of section 2.4. This condition was that V{{ki}) should stay invariant under certain 
transformations of the ki which kept the {V{ki)} constant. It allowed phase-space 
integrations over rapidity to be carried out analytically. 

Nevertheless we can still use the ‘general’ method for the elimination of sub¬ 
leading logs, as outlined in section 2.4. It requires that we be able to evaluate the 
observable accurately even when it takes values much smaller than the smallest ma¬ 
chine representable (double precision) floating point number, so that we can take the 
limit ttg 0 with constant agL (which implies L 00 ). 

It turns out that to succeed in doing this (without recourse to arbitrary precision 
arithmetic) we must carry out an analytical analysis of the observable’s dependence 
on multiple soft and collinear emissions. One might raise the criticism that this 
negates the philosophy of our numerical approach. Such a criticism would be only 
partially justihed since in a traditional resummation approach once one has deter¬ 
mined this dependence, one still has to understand what transformations are needed 
in order to write it in a factorised form (c/. the discussion in the introduction), and 
then one has to carry out the inverse transformations, both operations involving 
considerable work. 

4.3.1 Soft and collinear analysis 

It is useful to start by recalling the dehnition of the Durham jet hnding algorithm [26], 
given in terms of a resolution parameter |/cut as follows: 

1. For all pairs of (pseudo)particles i,j calculate 

E]){1 - cos%) 

= - Q2 - ■ (4-3) 

2. If all i/ij > I/cut stop. The number of jets is then dehned to be equal to the 
number of (pseudo)particles left. 

3. Otherwise recombine the pair with the smallest i/ij into a single pseudoparticle 
(there are various schemes according to whether one wants the recombined 
particle to be massless, see for instance [27]). 
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4. Go back to step 1. 


The three-jet resolution parameter, 1 / 3 , is the maximum value of |/cut that leads to a 
3-jet event. 

To understand the soft and collinear limit of our observable, we shall suppose 
that we have an event with hard partons a and b and soft-collinear partons 1 ... n. 
To start with we should examine the value of the for different kinds of pairings 
of particles: 

• If i and a are in the same hemisphere i/ia — Ef9f ~ 

• If i and j are in different hemispheres then yij ~ 4i?|, where = min(i?j, Ej). 

• If i and j are in the same hemisphere then ~ E‘^\ 6 i — 6 j\‘^, where we dehne 
the angles as vectors in the transverse directions. 

In these relations we have dehned 6 i to be the angle between parton i and its nearest 
hard parton (the vector indicates the azimnthal direction), and ku to be the relative 
transverse momentum. The centre of mass energy Q is taken to be 1. 

To determine the valne of y^ for an arbitrary ensemble of soft and collinear 
particles it is benehcial to start hrst with some simple cases: if we have just a single 
emission then y^ = kf of that emission (hence the alternative name — kt algorithm). 

Now let us consider two emissions, labelled 1 and 2. We take ku > kt 2 , the 
non-trivial case being when kt 2 and ku are however of the same order. We shall 
assnme that the two emissions are widely separated in rapidity, as is appropriate 
when considering only NLL terms. There are then four relevant conhgurations, 
illustrated in hgure 6 : 

i) If 1 and 2 are in opposite hemispheres (as in hgnre) then yi 2 = 4 min(i?f, i?|) 3> 
Via > y 2 b, so parton 2 gets recombined with b. Since now the only soft-collinear 
parton left is 1 we get y^ = kf^. 

a) If 1 and 2 are in the same hemisphere (say the hemisphere of a) and di ^ 62 
(which implies Ei i? 2 ), then we get yi 2 — yia = kfi, while y 2 a = k ^2 which 
is smaller. So once again parton 2 recombines with a hard parton and we are 
left with just parton 1 , which then leads us to y^ = k^^. 

Hi) If 1 and 2 are in the same hemisphere (again the hemisphere of a), bnt 61 62 

(which implies Ei 3> E 2 ), then we have yi 2 — y 2 a = k‘f 2 , while yia = k^ which 
is larger. So the recombination that will take place is either 12 or a2. Which 
one depends on whether 1 6*2 — 6^1 p > 6 ^ 2 . If 6*1 ■ 6^2 < 0 then the inequality is 
satished, and we have a a 2 recombination, so that just as before y^, = k^^. 
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Figure 6: Representation of the clustering in the Durham algorithm, for various config¬ 
urations. Both gluons are soft and collinear, with kti > kt2 but of the same order. When 
they are in the same hemisphere, their angles are assumed to be strongly ordered (though 
for clarity only weak ordering is represented in the figure). The hoop indicates which pair 
of particles clusters first. 

iv) If 1 and 2 are in the same hemisphere and 6 *i -C 6 * 2 , as in zn), bnt with 01 - 6*2 > 0, 
then yi 2 < y 2 a and (at last a non trivial resnlt!) we recombine partons 1 and 
2. This gives ns a psendoparticle with the energy of parton 1 (recall Ei^ E 2 ) 
and sqnared transverse momentnm = \kti + /Ci 2 p- The valne of r /3 is then 
just set by the of the pseudoparticle. We note that since 9i ■ 62 > 0, ktp is 
always larger than k^. 

We can see from this analysis that item 3 of the conditions for applicability in ap¬ 
pendix A is violated: exchanging the rapidities of the two gluons, while keeping their 
transverse momenta and azimuths constant (he. exchanging conhgurations ii and 
iv), leads to a change in the value of 7 / 3 . 

The above two-gluon analysis also points the way to a general algorithm suitable 
for taking the limit of very small values of y^. In order to be able to represent all 
quantities on a computer in standard double precision we work in terms of rapidities 
r]i and rescaled transverse momenta Ki = ku/ maxj{fcij}. We further assume (as is 
appropriate at NLL order) that particles are all widely separated in rapidity from 
one another. The j /3 value for the Durham jet algorithm can then be determined as 
follows: 

1. Find the index I of the smallest of the {Ki}. 

2. Considering only the soft partons j in the same hemisphere as I and which 
satisfy kj ■ Kj >0, hnd the index J of the one with the smallest positive value 
of \rij \ — \rii\. If there are no soft partons j with both kj ■ Kj > 0 and \rij\ > \rii\, 
then let J = a or b, according to which is the collinear hard parton. 
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3. If J = a or J = 6 then just throw away parton /. Otherwise recombine / and J 
as follows: the rapidity rjp of the pseudoparticle is just set to rjj (this involves 
a mistake on rjp by an additive amount of order 1 — but that is subleading), 
while its transverse momentum is the vector sum, Kp = ki + kj. At NLL order, 
this procedure is appropriate for all standard recombination schemes. 

4. If only one soft pseudoparticle remains, P, then y^/ maxj{/c^j} = Kp. Otherwise 
go back to step 1 . 


4.3.2 Numerical results 

Using the above algorithm for calculating at arbitrarily small values we can now 
numerically determine the function The results are shown in hgure 3 to¬ 

gether with those for the thrust major and the oblateness.® 

For phenomenology one also needs to know the resummation for the simple 
observable related to y^. Since y^s is just the (normalised) squared transverse mo¬ 
mentum, |/ 3 s = /Cj/Q^, we have from (2.14) and (2.15) that 


Ss(|/3) 


— Q-Rsiys) 


RsiUs) — Cp 



dk^ as{kt) 
kt vr 



(4.5) 


This corresponds to what was calculated in [14], which contains the explicit NLL 
form for the result of the integration. The full NLL resummed prediction for the 
Durham 1/3 can then be obtained from (2.8), with Ry^{R') as shown in figure 3. 

For completeness, in hgure 7 we show the hnal perturbative result for the NLL 
resummed Durham y^ distribution, matched to the NLO hxed order results. We 
choose here as a recombination scheme the PO scheme (see for instance [27]), in 
which the energy of the recombined pseudo-particle equals the sum of the energies 
of the two (pseudo)-particles, while the three-momentum is rescaled so as to keep 
pseudo-particles massless. Choosing a different recombination scheme will affect only 
the large y^, region and subleading logs. In the same hgure we also show the matched 
resummed curve for the Cambridge algorithm [28], which is also in the ‘family’ of 
kt algorithms — however it has a diherent clustering sequence which leads to the 
property that for emissions widely separated in rapidity, the 1/3 value is simply the 
maximum of the emitted squared transverse momenta.^ As a consequence R{R') = 1, 
so that the full resummed distribution is simply given by S^. 

®We note that the above analysis also makes it quite simple to determine analytically: 


^ 1 d(j) 7^" 1 kfi+k^2 + ‘^ktikt2COS(j> tt^ 


(4.4) 


where the leading factor of 1/4 comes because only a quarter of the time are the particles in the 
same hemisphere with 9i <C 02 - 

^We would like to thank Yu.L. Dokshitzer and B.R. Webber for bringing this to our attention. 
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Figure 7: The full SL-resummed distributions for the Durham and Cambridge three-jet 
resolution including second order matching (log i?-scheme). 


5. Conclusions 

In this paper we have presented a general method for the numerical calculation of 
the only non-trivial class of next-to-leading logarithms in (global) event shapes. For 
most event shapes it means that the calculation of the resummed distribution (in the 
n-jet limit for an observable which has non-zero values only for ensembles of n -|- 1 
or more particles) can be reduced to the following three tasks: 

• The evaluation of the dependence of the observable on a single soft and collinear 
emission. 

• The calculation of a resummed distribution, based on that single-emission 
dependence, through eqs. (2.14) and (2.15) (or their extensions in situations 
with more than two hard ‘legs’, and/or incoming hadrons). 

• The provision of a computer subroutine which calculates the value of the ob¬ 
servable for an arbitrary set of four-momenta. One can then use the algorithm 
of this paper to determine a single logarithmic function T{R') which accounts 
for the observable’s non-trivial dependence on multiple emissions, and multi¬ 
plies Tjg to give the resummed distribution to full NLL accuracy. 

The method thus provides both a very simple, but also general way of obtaining NLL 
resummed distributions for hnal-state observables. 
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As a cross check, it has been tested against a range of previously studied ob¬ 
servables (including some in processes other than which for brevity we have 

not shown here, but which are listed at the end of appendix B), and consistently 
reproduces the known analytical results. 

It has also been applied to observables that up to now had proved beyond the 
scope of existing analytical methods, allowing us to present the first fully NLL re¬ 
summed distributions for the thrust major, the oblateness and the three-jet resolution 
1/3 in the Durham algorithm. For the three-jet resolution some analytical analysis 
was required in order to provide a subroutine capable of giving an accurate evalua¬ 
tion of 1/3 even in the limit of very small y^. That analysis also has some intrinsic 
interest for the light it casts on the functioning of the Durham jet algorithm. 

It should be emphasised that the philosophy of the approach presented here dif¬ 
fers entirely from that of a Monte Carlo event generator (such as Herwig [29] or 
Pythia [30]). An event generator may well reproduce many of the leading and sub¬ 
leading logarithms for event shape distributions, however even at parton level there 
will also be contamination from potentially spurious NNLL and non-perturbative 
terms, and it is not currently possible to match to exact next-to-leading order calcu¬ 
lations. 

On the other hand the results produced by our algorithm are indistinguishable 
from those of fully analytical methods, allowing the same procedures of truncation 
at pure LL and NLL terms, and straightforward matching to fixed leading and NLO 
calculations. As such, our method is an essential development should one wish at 
some stage to fully automate the calculation of resummed distributions, in a manner 
analogous to what has become standard for hxed-order predictions of non-inclusive 
observables. 
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A. Conditions for applicability 

Various steps in this article rely on the observable under consideration being ‘suit¬ 
able’. Here we attempt to dehne what is meant by this. The hrst two conditions 
are required for the exponentiation of the double logarithms. The third condition is 
specihc to the second (event-shape specihc) method for eliminating subleading effects 
(discussed in section 2.4): 

1. If one is resumming from a Born level consisting of n coloured hard partons 
(incoming or outgoing), then for any conhguration consisting of just n hard 
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partons the observable must have a constant value Vq independent of the hard 
conhguration. Furthermore for a general conhguration of n + 1 hard partons 
the observable must have a value different from Vq. Without loss of generality, 
one can redehne the observable such that Vq = 0 . 

So for example, we can consider the thrust in the 2-jet limit, or the thrust 
minor in the 3-jet limit, but not the thrust-minor in the 2-jet limit, nor the 
thrust in the 3-jet limit. 

2. If we select m soft, collinear momenta (with respect to any of the n hard 

partons) such that V{ki) ~ V{k 2 ) V{km), then the observable must 

satisfy the condition V{ki ,... , km) ~ V{ki). This is necessary in order for the 
double logarithms to exponentiate, and excludes for example observables based 
on the JADE jet-clustering algorithm [31]. 

3. Given m soft and collinear (SC) partons, the value of the observable V{ki ,... , km) 
must remain unchanged (to within corrections of relative order of the softness 
or collinearity) if we vary the rapidity and absolute transverse momentum (but 
not the azimuth) of any of the m SC partons subject to the following restric¬ 
tions: 

• the C(A;j) are all kept constant; 

• each SC parton remains collinear to the hard parton to which it was 
originally collinear. 

This condition is necessary if one is to use method 2 of section 2.4 for the 
elimination of subleading logs. We note that in general it is satished by event 
shapes, but not for example by the jet rate with the Durham algorithm. 

In addition to the above conditions, we note that for the simple observable to have a 
straightforward resummation in terms of independent emissions (c/. eqs. (2.12), (2.14) 
and (2.15)) it is necessary that it be global [12]. 


B. Definition of observables 


For completeness we recall here the dehnition of all the event shape variables that 
we refer to in the main text (except for whose dehnition is explicitly needed in 
section 4.3 and so is directly given there). 


• Thrust T, Thrust major Tm, Thrust minor T^, Oblateness O 


The thrust is dehned as 


T 


Ei \Pi 


max 




\Pi 


■ n\ 


(B.l) 
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where pi are the three-momenta of the outgoing particles. The thrust axis is 
the unit vector ut which maximises the sum in (B.l) and gives therefore the 
direction along which the projection of momenta is maximal. 

The thrust major is then dehned similarly, but the maximisation procedure is 
restricted to unit vectors perpendicular to Ht 


Tm = 


max 


E. \p^\ -->=0 


\Pi ■ n 


(B.2) 


The vector htm which maximises the sum in (B.2) dehnes the thrust major 
axis. 

Given the thrust and the thrust major axes, the thrust minor is dehned as 

Tm = ^^1^1 V \Pi -nTj, ut^ = nTX , (B.3) 

E* m V 

and the oblateness as 


O — Tm — Tm . 


(B.4) 


• The broadenings: 

The broadenings measure the transverse size of the jets. The plane perpendic¬ 
ular to the thrust axis through the origin divides an event in two hemispheres, 
which are usually called left (L) and right (R) hemisphere. One dehnes hrst 
the right {Bfj) and left broadening as 


Br/l 


1 

2 Ei \Pi 


’ 

ieRIL 


(B.5) 


with pti = Pi X Ut- 

One can then dehne the total (Br) and wide {Bw) jet broadening as 


Bt — Br + Br , 

Bw = max{R/j, Br} ■ 


(B.6) 


As well as having been applied to the observables dehned above, our method 
has also been successfully tested against known analytical results for several other 
variables: the D-parameter in 3-jet e’''e“ events [9], the current-jet broadening (with 
respect to the photon axis) in 1 -|- 1-jet DIS events [6], the out-of-plane momentum 
in 2 -|- 1-jet DIS events [11] and the out-of-plane momentum for W/Z -|-jet events in 
hadron-hadron collisions [10]. 
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C. Analytical ingredients 

Here we summarise all the analytical ingredients needed for the calculation of the 
resummed distributions of the three observables studied in section 4. This means 
NLL expressions for the simple radiators and derivatives (the formulae are all well 
known in the literature), as well as the hxed order expansions. 

Computer subroutines which interpolate tabulated numerical values for the func¬ 
tions are available on request from the authors, as are example programs for 

calculating the complete resummed, matched distributions. 


C.l Thrust major and oblateness 

For the thrust major and oblateness, the simple radiator, Rg (dehned in (4.2)) is 
given to NLL order by Rs{v) = Lgsi{asL) + gs 2 {asL) with 

Cf 


gsiiasL) = 


(-2A-ln(l-2A)), 


and 




2A) + 


+ 


7rf3o\ 

CpK 2A + (1 -2A)ln(l -2A) 


(C.l) 


CfPi 

TTpi 


1 - 2A 

2A + ln(l - 2A) 


-ln2(l 
2 ^ 


2A) , (C.2) 


1 - 2A 

where A = asj3oL and L = In 1/u. The coefficients of the /5-function are given by 
lie A -2nf 17CI - hCAUf - SCpUf 


Po — 


Pi = 


127r ’ 247r2 

and the constant relating the gluon Bremsstrahlung scheme [23] to the ms is 

Q7 tt^A 5 
18~y)~9^^' 

Also useful is the formula for R' (at NLL equal to /?(): 

4A 


K = Ca 




2A 


(C.3) 


(C.4) 


(C.5) 


For matching to fixed order calculations, it is necessary to have the coefficients 
of the hxed order expansion of the full resummed distribution. They are given in 
table 1, using the notation of [1]. 

Tables 2a and 2b show the two subleading coefficients G 21 and C 2 for the thrust 
major and the oblateness respectively. They are obtained from the hxed order inte¬ 
grated distribution (calculated with EVENT2 [16]) after subtracting the leading and 
next-to-leading logs. The hnal errors shown are obtained by adding in quadrature 
the mean error and the maximal discrepancy between the mean coefficient and the 
coefficient obtained in a single ht. The errors are to be considered as indicative only, 
because of the difficulty in estimating the systematic uncertainty associated with the 
choice of ht range. 
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Gi 2 — —ACp 

Gii = Cf {6-8\n2) 

Cl = C'ir(61n2-41n2 2 + 7r2-f) 

G23 = -CpifCA-fTnUf) 

G 22 = C'^(64^2C'i.+ (|7r2-f-fin 2)^^+(I + fin 2) T^nj) 

Table 1: Logarithmic coefficients Gnm and the constant Ci for thrust major and oblate¬ 
ness; for Tm one has T 2 = —0.46851 ± 0.00017, for the oblateness T 2 = 0.37799 ± 0.00029. 


a) Thrust Ma^ 

or 

Fit-range 

021 

02 

[-9;-5] 

-20.3 ±2.5 

154.0 ± 14.8 

■ 1 ' 
cc 

1 

-20.7 ±3.3 

156.2 ± 19.3 

[-9;- 6 ] 

-21.0 ±2.9 

160.4 ±19.2 


-20.6 ±2.9 

156.9 ±18.1 


b) Oblateness 


Fit-range 

021 

O 2 

1 - 8 ; -4] 

75.3 ± 1.0 

166.1 ±4.6 

' 1 ' 
-u 

1 

75.6 ± 1.1 

164.9 ±5.3 

■ 1 ‘ 
00 

1 

70.8 ±2.9 

191.6 ± 16.5 


73.9 ±3.5 

174.2 ± 19.5 


Table 2: Fits of G 21 and C 2 for the thrust major and the oblateness. 


C.2 Three-jet resolution 


For the Durham (and Cambridge) three jet resolution parameter, the simple radiator, 
Rg, defined in (4.5), is given to NLL order by Rgiy) = Rgsi^a^iL) + gs 2 {aiiL) with 


[-^ - In (1 - A)] , 


(C. 6 ) 


^ + 27r2/32(l - A) 


+ 


CfPi 

L 


A -I- ln(l — A) 1^2 


1-A 


--ln^(l-A) 


and, as before, A = agPoL- R' is given by 

Cf a 


R' = 


(C.7) 


TT/do 1 — A 

The coefficients of the fixed order expansion of the resummation are given in table 3. 
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Table 3: Logarithmic coefficients Gnm and the constant Ci for the three-jet resolution 
parameter. In the Durham algorithm = —'7r^/32, while in the Cambridge algorithm 
T 2 = 0 . 

a) Durham b) Cambridge 


Fit-range 

^21 

^2 

[-13;-8] 
[-13;-9] 
[-12;-8] 

1.9 ±0.1 

2.3 ±0.2 

1.9 ±0.2 

-30.5 ± 1.2 

-34.1 ± 1.6 

-30.3 ± 1.4 


2.1 ±0.3 

-31.6 ±2.9 


Fit-range 

^21 

^2 

[-13;-8] 
[-13;-9] 
[-12;-8] 

-7.4 ±0.1 

-7.0 ±0.1 

-7.4 ±0.1 

19.3 ± 1.0 

15.8 ±1.3 

19.4 ±1.1 


-7.2 ±0.3 

18.2 ±2.6 


Table 4: The coefficients G 21 and C 2 for the three-jet resolution in the Durham and 
Cambridge algorithms. 


From fixed order results we also extract the subleading coefficients, shown in 
tables 4a and 4b for the Durham and Cambridge algorithms respectively. We follow 
the same £t procedure as for the thrust major and the oblateness. Note that these 
coefficients are obtained in the F^O-scheme (see for example [27]), and that choosing 
a different recombination scheme will affect their value. 
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